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ABSTRACT 

We build a theoretical model to study the origin of the globular cluster metallicity bimodality in the hierar- 
chical galaxy assembly scenario. The model is based on empirical relations such as the galaxy mass-metallicity 
relation [O/H] -M stal - as a function of redshift, and on the observed galaxy stellar mass function up to redshift 
z ~ 4. We make use of the theoretical merger rates as a function of mass and redshift from the Millennium 
simulation to build galaxy merger trees. We derive a new galaxy [Fe/H] - M stal - relation as a function of red- 
shift, and by assuming that globular clusters share the metallicity of their original parent galaxy at the time 
of their formation, we populate the merger tree with globular clusters. We perform a series of Monte-Carlo 
simulations of the galaxy hierarchical assembly, and study the properties of the final globular cluster popula- 
tion as a function of galaxy mass, assembly and star formation history, and under different assumptions for the 
evolution of the galaxy mass-metallicity relation. The main results and predictions of the model are the follow- 
ing. 1) The hierarchical clustering scenario naturally predicts a metallicity bimodality in the galaxy globular 
cluster population, where the metal-rich subpopulation is composed of globular clusters formed in the galaxy 
main progenitor around redshift z ~ 2, and the metal-poor subpopulation is composed of clusters accreted from 
satellites, and formed at redshifts z ~ 3 - 4. 2) The model reproduces the observed relations by Peng et al. 
(2006) for the metallicities of the metal-rich and metal-poor globular cluster subpopulations as a function of 
galaxy mass; the positions of the metal-poor and metal-rich peaks depend exclusively on the evolution of the 
galaxy mass-metallicity relation and the [O/Fe], both of which can be constrained by this method. In partic- 
ular, we find that the galaxy [O/Fe] evolves linearly with redshift from a value of ~ 0.5 at redshift z ~ 4 to 
a value of ~ 0.1 at z = 0. 3) For a given galaxy mass, the relative strenght of the metal-rich and metal-poor 
peaks depends exclusively on the galaxy assembly and star formation history, where galaxies living in denser 
environments and/or early types galaxies show a larger fraction of metal-poor clusters, while galaxies with a 
sparse merger history and/or late type galaxies are dominated by metal-rich clusters. 4) The globular cluster 
metallicity bimodality disappears for galaxy masses around and below M stal - ~ 10 9 Mq, and for redshifts z > 2. 
Subject headings: Galaxies: star clusters: general - Galaxies: formation - Galaxies: evolution- Galaxies: stellar 
content - Galaxies: structure - Galaxy: globular clusters: general 
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1. INTRODUCTION 

Globular cluster (GC) systems in galaxies have become 
a useful tool to study the mechanisms of galaxy formation. 
Thanks to a rise in the level of details in observations, now 
we can gain insight into the colour, metallicity and abun- 
dance gradients of such systems for a large number of galax- 
ies, and build statistically solid scaling relations between GC 
and galaxy properties. 

GCs are for the most part old objects, with ages estimated 
to be > 10 Gyr (Brodie et al. 2005, Strader et al. 2005, Peng 
et al. 2006). Therefore, not only they have survived any vio- 
lent event in the assembly of their host galaxy, but they also 
provide a chemical record of the galaxies where they were 
formed (Pota et al. 2012). Coupled with the fact that they 
are very luminous, they make for excellent probes of the fos- 
sil records of galaxies and shed light on the mechanisms of 
galaxy assembly and star formation history. 

Of particular interest is the metallicity distribution of GCs 
in galaxies. Galaxies of all morphologies have a GC popula- 
tion with an average metallicity that correlates with the galaxy 
stellar mass or luminosity (as first shown by Brodie & Huchra 
1991; Lotz et al. 2004, Peng et al. 2006). In addition, most 
galaxy GC systems exhibit a colour bimodality (Zepf & Ash- 
man 1993, Ostrov et al. 1993, Whitmore et al. 1995, Elson & 
Santiago 1996, Peng et al. 2006, Spitler et al. 2006, Strader 



et al. 2006, Larsen et al. 2001). This is driven by a metal- 
licity bimodality, with bluer GCs being more metal-poor and 
redder GCs being more metal-rich, while both populations are 
old (> 10 Gyr) (Forbes et al. 2001, 1997a,b, 201 1, Peng et al. 
2006, Strader et al. 2005, 2006, Cote et al. 1998, Puzia et al. 
2005, Pierce et al.2006, Brodie et al. 2005, Brodie & Strader 
2006). Although there is some debate in the literature (see 
Yoon et al. 2006), this result has been confirmed spectroscop- 
ically (Brodie et al. 2005, 2012, Cohen et al. 2003, Strader 
et al. 2005, Alves-Brito et al. 2011, Usher et al. 2012). The 
well-defined metal-rich and metal-poor GC sequences sepa- 
rately follow two galaxy stellar mass - GC metallicity rela- 
tions [Fe/H]cc - M stal , of which the metal-rich one is stronger 
and tighter, while the metal-poor one is weaker and exhibits a 
larger scatter (Larsen et al. 2001, Strader et al. 2006, Peng et 
al. 2006, Cote et al. 1998). There is increasing evidence that 
these features of GC systems are universal, from giant ellip- 
ticals to dwarfs (Strader et al. 2006), although some galaxies 
show an even more complex situation, with multiple metal- 
licity peaks (see for instance Peng et al. 2006, Blom et al. 
2012). 

Intriguingly, the differences between blue/metal-poor and 
red/metal-rich GCs also extend to their dynamical properties, 
as shown in recent observations (Pota et al. 2012). The two 
subpopulations have different spatial distributions inside the 
host galaxies: the metal-rich GCs are more centrally concen- 
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trared, with a radial distribution profile that follows closely 
the spheroidal stellar component of the galaxy, while the 
metal-poor GCs show a more extended distribution, and is 
likely physically associated with the stellar halo (Bassino et 
al. 2006, Goudfrooij et al. 2007, Peng et al. 2008, Forbes 
et al. 2012, Pota et al. 2012; the M87 data of Strader et 
al. 2011 represent the best example of the close spatial cou- 
pling of metal-rich GCs with galaxy starlight and the more 
extended distribution of metal-poor GCs). Correspondingly, 
the kinematics of the metal-rich subpopulation follows that 
of the main stellar component, including rotation (Strader et 
al. 2011), while the metal -poor subpopulation shows larger 
velocity dispersion and small or null net rotation. 

A scenario has been proposed where GCs are formed in 
gas-rich (major) merger events; at high redshift (z > 4-5), 
early mergers of smaller hosts produce metal-poor GCs, while 
later mergers of more evolved galaxies in high density envi- 
ronments produce metal-rich GCs (Muratov & Gnedin 2010, 
Kravtsov & Gnedin 2005, Bekki et al. 2007, 2008). These 
models however encounter a number of problems; there is no 
clear prediction about any metallicity bimodality or galaxy 
mass-GC metallicity relations, and the resulting ages of the 
metal-rich GCs are too young (Muratov & Gnedin 2010), an 
ad hoc mechanism is needed to shut off blue/metal-poor GC 
formation (Bekki et al. 2008, Beasley et al. 2002), and an 
analysis of the observed GC abundance and metallicity gra- 
dients is not compatible this kind of formation mechanism 
(Arnold et al. 2011). 

Alternatively it has been proposed that, rather than origi- 
nating from two main epochs or modes of GC formation, the 
GC chemo-dynamical bimodality can stem from the galaxy 
assembly history, without invoking mergers as the GC for- 
mation mechanism. In this scenario the metal-rich GC sub- 
population is formed together with the bulk of the galaxy 
stellar component in an early violent dissipative phase, and 
during a later slower phase the metal-poor GC subpopulation 
is accreted, via minor mergers (Forbes et al. 2011, 1997a,b, 
Arnold et al. 201 1, Masters & Ashman 2010), or via stripping 
of GCs from satellites (Cote et al. 1998, 2000). The main 
difference with the merger scenario is that GCs of different 
metallicities are formed in different galaxies, and then brought 
together by galaxy assembly, rather than being formed in the 
same galaxy at different stages of the galaxy evolution. In this 
work we call this the "assembly scenario". 

1.1. This work 

In this work we want to put the "assemby scenario" in the 
context of the hierarchical structure formation theory, and in- 
vestigate whether the GC metallicity bimodality indeed orig- 
inates from the hierarchical nature of galaxy assembly. In 
other words, is the GC metallicity bimodality a natural pre- 
diction of hierarchical clustering? 

To answer this question, we build a model to produce the as- 
sembly history of galaxies and their GC population, in a series 
of Monte-Carlo simulations. We base our model galaxy prop- 
erties on observed scaling relations as a function of redshift, 
such as the galaxy mass-metallicity relation [O/H] -M star re- 
lation, and the galaxy stellar mass function. We assume that 
galaxies at z = were formed through a combination of local 
(in-situ) star formation and accretion of satellite galaxies in a 
series of merger episodes spanning the lifetime of the galaxy; 
the merger rates are obtained from the Millennium simula- 
tion. We populate each galaxy in the merger tree with GCs, 
assuming that they share the metallicity of the main stellar 



component of their parent galaxy at the epoch of their forma- 
tion. When a satellite is accreted, so is its GC population. 

We investigate under what conditions the final GC popu- 
lation shows the metallicity bimodality, and follows the ob- 
served metal-rich and metal-poor galaxy stellar mass - GC 
metallicity relations, as well as the observed galaxy mass - GC 
number abundance relation (Peng et al. 2006, 2008, Strader 
et al. 2006). 

The novelty of this analysis is that it provides constraints 
and predictions 1) on the galaxy [Fe/H] -M star relation as a 
function of redshift, 2) on the galaxy assembly and star forma- 
tion history, and 3) on the evolution of the GC bimodality, and 
ultimately it presents a method to test the hierarchical galaxy 
formation. 

This paper is organised as follows. In Sections 2 and 3 we 
present the model: in Section 2 we describe the galaxy assem- 
bly and the globular cluster formation, and the Monte Carlo 
simulation; in Section 3 we present the derivation of the fidu- 
cial galaxy mass-metallicity relation. In Section 4 we present 
our results for the globular cluster metallicity distribution and 
its implications to constrain galaxy formation, and in Section 
5 we discuss our findings. Section 6 is a summary of our con- 
clusions. 

2. THE MODEL: GALAXY ASSEMBLY AND 
GLOBULAR CLUSTER FORMATION 

Consider a galaxy of stellar mass Mo at redshift z = 0. This 
object represents the final stage of a merger tree, i.e. a sys- 
tem of independent progenitor galaxies which were accreted 
and contributed to all the mass components (dark matter, gas, 
stars, GCs) that now characterise the galaxy. At any given 
time, we identify the main progenitor in the merger tree as 
the most massive galaxy that is present in the tree, while we 
(improperly) call satellites the rest of the objects. For any 
given galaxy at z = 0, we build a Monte Carlo simulation 
with N realisations of the merger tree, i.e. N different assem- 
bly histories. We performed numerical tests on N in a range 
N = [10, 10 6 ], finding convergence of our results for N > 100. 
The plots in this work are made with N = 10 3 . 

Galaxies evolve depending on the mass of the host dark 
matter halo and on the density of the surrounding environ- 
ment. In the hierarchical clustering framework smaller ob- 
jects virialise earlier (see for instance Frenk & White 2012), 
so they contain older, metal-poor stars. Their cycle of star for- 
mation and feedback is less efficient, and supernovae winds 
are more effective in expelling metals from the galaxy, factors 
that contribute to slow down the rise of metallicity in their 
stellar populations. At the same time, in more massive galax- 
ies the deeper potential wells render supernovae winds less 
effective in expelling metals, and the enhanced ability to re- 
tain gas allows for sustained star formation and more stellar 
generations. As a consequence, at all redshifts a monothonic 
positive mass-metallicity correlation M star -[Fe/H] is in place 
for all galaxies in the merger tree. The derivation of this rela- 
tion from the observed [0/H]-M star relation will be described 
in detail in the next Section. 

We assume the appearance of a globular cluster popula- 
tion in a galaxy is an event of a relatively short duration, 
and in general not associated with the quiescent star forma- 
tion phase, but indicative of a particularly intense evolution- 
ary phase. This assumption is sustained by a number of ob- 
servations. First, the observed masses of GCs can reach up 
to M ~ 10 6 M , requiring very intense bursts of star forma- 
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tion. Secondly GCs, which in general are well described by 
single stellar populations (SSP, i.e. coheval ensembles of stars 
that share the same metallicity), are for the most part old ob- 
jects, with ages > 10 Gyr (as referenced in the Introduction). 
Fittingly, the observed GC ages put the epoch of their appear- 
ance squarely at the peak of the cosmic star formation history, 
determined to be at redshifts z > 2 - 4 (Hopkins & Beacom 
2006, Bouwens et al. 2009). Third, the mean metallicity of 
the GC population is observed to be higher in more massive 
galaxies, with a M star - [Fe/H]cc parallel to that of the galaxy 
mass-metallicity relation (see for instance Larsen et al. 2001, 
Peng et al. 2006, Cote et al. 1998), in support of the idea 
that the GC population is closely related to the main stellar 
component and is similarly affected by halo mass and envi- 
ronment, i.e. more massive galaxies form their bulk of their 
stars later, from more enriched gas. 

We assume that the old globular cluster population (ages 
> 10 Gyr) were formed in galaxies at the peak of their star 
formation activity, likely L* galaxies, at all redshifts z > 2. 
Note that, following this assumption, the older globular clus- 
ters formed at higher redshifts in smaller systems, and with a 
lower metallicity. This is in accord with observations, that es- 
timate metal-poor GCs to be about 1 - 2 Gyr older than metal- 
rich GCs (Dotter et al. 2011, Puzia et al. 2005, Woodley 
et al. 2010), although the precision of the age measurement 
for extragalactic GCs is too low to confidently discriminate 
ages differences at this level (Strader et al. 2005). We assume 
as likely candidates for the formation sites of GCs either the 
massive star-forming clumps observed in high redshift (z > 2) 
galaxies (Shapiro et al. 2010), or the central regions of galax- 
ies subject to episodes of violent dissipative collapse. In both 
cases, the dumpiness and turbulence of the gas plays a funda- 
mental role in boosting the star formation and producing GCs, 
along with the galaxy main stellar component (Shapiro et al. 
2010). 

The frequency of globular clusters Tn is defined as the num- 
ber of GCs per unit galaxy mass of 10 9 M , and at z ~ it 
is constrained by observations (Peng et al. 2008; see also 
Spitler et al. 2008, Rhode et al. 2007, Rhode 2012). In the 
galaxies in the merger tree, Tn depends on the interplay of 
different factors, like the mean gas density (which depends 
on the depth of the galactic potential well), the metallicity, 
the feedback regime, and the competing 'regular' star for- 
mation that feeds the main stellar component. In lack of 
other observational constraints, we assume that the redshift 
z <~ observed relation T N -M star holds at all redshifts, 
so that the total number of local GCs that each galaxy pro- 
duces is Ngc = TN(M star ) x M star . In addition, galaxies below 
M star = 10 9 M© = M m i n do not form globular clusters, consis- 
tently with the observed T N -M star relation (Peng et al. 2008), 
which yields N GC < 10 for M star ~ 1O 9 M ( see also Muratov 
& Gnedin, 2010). We also assume that, once formed, the lo- 
cal GCs stabilise themselves in dynamical equilibrium with 
the galaxy, and therefore remain kinematically coupled with 
the main stellar component. 

We follow the evolution of the merger tree from redshift 
Z[ n , when the galaxy main progenitor forms its local globular 
clusters. We assume the epoch is z ln ~ 2. The main progeni- 
tor is likely to have a rich gas component and is near the peak 
of its star formation history; it has a stellar mass Mi and a 
mean total metallicity Zi, which follows the galaxy mass - 
metallicity relation M star - [Fe/H] at z ~ z ln . In each Monte- 
Carlo realisation, we assume that the metallicity of the locally 



formed GCs is peaked aroung Zi, with a gaussian distribu- 
tion with a = 0.2 (consistent with Cote et al. 1998, Bekki et 
al. 2008), which takes into account a non-instantaneous mix- 
ing of the metals, and the fact that the GC formation covers a 
short but finite time-span, in which the mean galaxy metallic- 
ity can vary. The number of local GCs in the main progenitor 
isN GC = T N (M 1 )xM 1 . 

The main progenitor is the most massive galaxy in the 
merger tree, and is the last one to have its globular cluster 
population in place. At this point in time, the satellites in 
the merger tree, which by definition have masses M; < Mi, 
have already formed their own GCs, and we put such epoch 
of formation around z ~ 3-4 (see also Shapiro et al. 2010). 
Each satellite metallicity follows the galaxy M star - [Fe/H] at 
that epoch, and in each satellite of mass Mj the GC metallic- 
ity is peaked around the current galaxy mean metallicity Zi, 
while the number of GCs is N G ci = T N (Mj) x Mj. We assume 
a gaussian distribution of the GC metallicities in each satel- 
lite, peaked around Z; with a = 0.3, which takes into account 
the combined effect of the scatter in the GC metallicity in each 
satellite (<~ 0.2), plus an additional uncertainty (~ 0.2) due to 
the scatter in the star formation histories of satellites (sensi- 
tive to environment for instance), which affect the satellite's 
metallicity and the exact epoch of GC formation. 

The main progenitor Mi evolves into the z = galaxy Mo 
through two main channels: by accreting stellar mass in the 
form of satellites, and by forming stars locally. If we define 
as Msf the mass in stars that are formed inside the main pro- 
genitor at any time after the GC formation (including merger- 
triggered star formation), then the stellar mass accreted from 
satellites is M sat = M -Mi -M SF . M sat is the sum of the stellar 
mass present in all satellites at redshift z ln , under the assump- 
tion that the satellite T N remains constant (i.e. satellites do not 
have a prolongued star formation history aftet the GC forma- 
tion). The ratios Mi /Mo and Msf/Mo are free parameters in 
the model, and they constrain the assembly and star-formation 
history. 

For each galaxy characterised by (Mo,Mi,Msf/Mi), we 
run a Monte-Carlo simulation of N realisations of the galaxy 
merger history, from z in to z = 0. In each run, we ran- 
domise the metallicity distribution of the main progenitor's 
GCs around Zi. We build the merger tree based on the 
observed stellar mass function (SMF) of Marchesini et al. 
(2009), and the theoretical merger rates obtained from the 
Millennium simulation (Springel et al. 2005, Fakhouri et al. 
2010). In each realisation, we randomise both the mass of 
the accreted satellites and the redshift of accretion; after sam- 
pling a random redshift in the interval [z ~ 4-0], we interpo- 
late the observed stellar mass function to that redshift, and we 
sample a random satellite mass from it with an acceptance- 
rejection algorithm. This provides us with a series of merger 
candidates; each of them is weighted with the mean merger 
rate, which represents the probability for a merger to happen, 
given the mass of the main progenitor M, the ratio between 
the masses of the satellite and the main progenitor e, and the 
redshift: 



dN m 
dedz 

where 



(M,e,z) = A 



=a(- 



M 



V1O 12 M 



exp 



(i+zY 



(i) 

the best-fit parameters are characterised 
as (a,/3, 7 ,7?) = (0.133,-1.995,0.263,0.0993) and 
(A,e ) = (0.0104,9.72 x 10" 3 ) (Fakhouri et al. 2010). 
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From Fig. (1) in Fakhouri et al. (2010) it is evident that the 
merger rate increases with increasing redshift and decreasing 
halo mass, and hugely favours small ratios e << 1 between 
the satellite's and the main progenitor's masses. 

At each timestep in our merger history, we add the weighted 
mass of each satellite to the stellar mass of the main progen- 
itor M, which grows in time, and we continue until the total 
accreted mass is equal to M sat . Each merged satellite carries 
a population of Ncci = Tn(M;) x M; globular clusters, with 
metallicity centered around Z; and randomised in each run. 
This completes one realisation in the Monte Carlo simulation 
and represent one of the N merger histories that we build for 
each galaxy. For each merger history, we obtain a total GC 
metallicity distribution that is given by the superposition of 
the contributions from the main progenitor and all the satel- 
lites. After N realisations, we produce a mean of the total GC 
metallicity distribution. 

In addition, we also explore a scenario where new globular 
clusters can be created in gas-rich merger events. In this case, 
we consider that at a redshift z new a merger event is charac- 
terised by a gas mass M eas that is turned into stars and globu- 
lar clusters, producing N new new globular clusters of random 
metallicity peaked around Z new (the metallicity of the gas), 
with a gaussian distribution with a = 0.2. 

For a given final galaxy mass, the model uses 2 free pa- 
rameters. The ratio Mi /Mo between the mass of the main 
progenitor at the epoch when it forms its local GCs and the 
final stellar mass of the galaxy is the assembly parameter, the 
ratio Msf/Mo between the mass of the stars formed locally in 
the evolving main progenitor after the GC formation and the 
final stellar mass of the galaxy is the star formation history 
parameter. 

3. THE MODEL: EVOLUTION OF THE GALAXY 
MAS S -METALLICITY RELATION AND GLOBULAR 
CLUSTER METALLICITY 

The main source of systematic uncertainty in the model 
comes from the redshift evolution of the galaxy mass- 
metallicity relation M stal -— [Fe/H]. Although this is in prin- 
ciple constrained by observations, we feel that there currently 
is a lack of consensus on the evolution of M stal - [Fe /H] at 
the level of precision required for this investigation. For this 
reason, to assign a metallicity to the galaxies in the merger 
tree and their globular cluster systems, we build a fiducial 
M star - [Fe/H] relation as a function of redshift, and we ex- 
plore the consequences of varying this relation on the model. 
As a sanity check, to obtain the metallicity of the satellite GCs 
we also make use of the total M star - [Fe/H]cc relation of Peng 
et al. (2006). 

The GC metallicity as a function of galaxy mass is provided 
in terms of [Fe/H], which is a proxy for the total metallicity 
(Fig. 14 of Peng et al. 2006, Shapiro et al. 2010). The galaxy 
metallicity on the other hand is often measured in terms of 
the quantity 12 + log(0/H); in particular we consider the rela- 
tions provided by Maiolino et al. (2008) up to z ~ 3.5 for the 
AMAZE (Assessing the Mass-Abundance redshift[-Z] Evolu- 
tion) program. In the lower redshift bin, this relation is con- 
sistent with the one provided by Tremonti et al. (2004) for a 
sample of 53000 galaxies in the Sloan Digital Sky Survey. 
The Maiolino relations can be parameterised as follows: 

12 + log(0/H) = -0.0864(logM stal -logM () ) 2 + ko , (2) 
and for redshifts z = (0.07,0.7,2.2,3.5) the parame- 
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FIG. 1. — The evolution of the galaxy mass-metallicity relation and the 
globular cluster metallicity relations from Peng et al. (2006). Thick lines: 
black average total GC relation, red: metal-rich GC relation, blue: metal-poor 
GC relation. Thin lines with shaded areas: the evolution of the galaxy mass- 
metallicity relation as obtained in this work (see text), with the la regions 
(z = 1 region omitted for clarity). The z = 0.1 dotted line is the relation we 
obtained from Tremonti et al. (2004). All the solid lines represent the relation 
we obtain from Maiolino et al. (2008). For the z = 3.5 relation we show 
the results for both Maraston (2005) and Brazual & Chariot (2003) stellar 
population models. 

ters are M = (11.18,11.57,12.38,12.76/12.87) and k = 
(9.04,9.04,8.99,8.79/8.9) (Maiolino et al. 2008). To ob- 
tain a [Fe/H] estimate from the quantity log(0/H) we need 
to establish the 12 + log(0/H) solar value, and the [O/Fe] or 
alternatively [a/Fe] values as a function of galaxy mass and 
redshift. These quantities are degenerate in producing the fi- 
nal [Fe/H]. 

The solar oxygen abundancy is determined to be 12 + 
log(0/H) = 8.66 (Erb et al. 2006), but other works put it 
at 12 + log(0/H) = 8.9 (as discussed for instance in Liu et 
al. 2008); unfortunately, the spread in the adopted value of 
the solar oxygen abundancy significantly increases the uncer- 
tainty in the calculation of the galaxy M star - [Fe/H] relation. 

In lack of direct spectroscopic observations, the determina- 
tion of [O/Fe] as a function of [O/H], of galaxy mass and 
of redshift depends on models of both stellar and galaxy evo- 
lution, and there is currently no consensus on the conversion 
[O/H] into [Fe/H] (A. Pipino et al. in preparation, and private 
communication). A determination of [a/Fe] as a function of 
galaxy mass at z ~ is provided by Thomas et al. (2005), for 
a sample of 124 early-type galaxies. The scatter is substantial, 
and the mass range does not include galaxies below 10 10 M Q . 
The relation is parameterised as follows: 

[a/Fe] = -0.459 + 0.062 logM star (3) 

If we use this prescription to convert the z ~ Maiolino and 
Tremonti relations, as [Fe/H] = [O/H] - [a/Fe], we obtain 
M star - [Fe/H] relations that are consistent (inside the scatter) 
with the one provided by Thomas et al. (2005) in the range 
M star = [10 10 - 1O 12 ]M , provided that [Fe/H] and [O/Fe] are 
reasonable proxies for [Z/H] and [a/Fe] respectively. The 
comparison yields values [a/Fe] ~ 0.1 for M stal - ~ 10 10 M Q 
and [a/Fe] ~ 0.18 for M star ~ 10 n M Q . 

We use the Maiolino et al. (2008) relations to obtain the 
[Fe/H]-M star relations at higher redshifts, but we need to 
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make an assumption about the redshift dependence of [a/Fe]. 
Such dependence is very uncertain and not all factors respon- 
sible for the variation of [a/Fe] are currently understood; for 
instance, a progressively top-heavy IMF at higher redshift 
would cause an excess of oxygen that would speed up the 
[a/Fe] evolution. For this reason, we choose to calibrate our 
[a/Fe] vs redshift relation empirically. Shapiro et al. (2010) 
use an estimated [a/Fe] ~ 0.3 to obtain a relation at redshift 
z - 2 from data of 12 + log(0/H) from Erb et al. (2006). A 
comparison with the z ~ 2 relation we obtain from Maiolino 
et al. (2008) via the Thomas et al. (2005) prescription, shows 
us that we need to assume that [a/Fe] evolves by 0.2 dex 
in order for the two relations to match. We then extrapo- 
late this evolution linearly with redshift, and obtain values 
[a/Fe] ~ (0.1,0.2,0.3,0.5) for redshifts z - (0, 1,2,3.5). 

We provide a rough estimate of the error in the [Fe /H] - 
M star relations from the scatter in the M stal - - 12+log(0/H) re- 
lation (~ 0.2 dex, Tremonti et al. 2004), the scatter in [a/Fe] 
at redshift (~ 0.1 dex, Thomas et al. 2005) and the uncer- 
tainty in the solar value of 12 + log(0/H) (~ 0.2 dex); we ob- 
tain an uncertainty a ~ 0.3 dex on [Fe/H] for any given stellar 
mass. Note that this estimate does not take into account the er- 
rors in the galaxy mass estimates, nor the error increase in the 
metallicity measurements at higher redshifts, and nonetheless 
a ~ 0.3 is of the same order of the systematic errors induced 
by our choice of the [a/Fe] evolution. We will explore the 
consequences of varying these relations in the next Section. 

In Fig. (Q]i we plot our fiducial galaxy [Fe/H]-M star re- 
lations up to redshift z ~ 3.5. On the same Figure, we plot 
the observed relations between the galaxy stellar mass and 
the GC metallicity [Fe/H]Gc _ M star obtained by Peng et al. 
(2006, their Fig. 14). These are shown as the straigh lines: 
black for the average GC metallicity in each galaxy, red for 
the metal-rich GCs, and blue for the metal-poor GCs. 

Once we have the galaxy [Fe/H] - M stal - relation in place as 
a function of redshift, we use it to assign a metallicity to all 
the globular clusters in the merger tree: the GCs formed in a 
galaxy of mass M stal - at a redshift z have a mean metallicity 
corresponding to the galaxy [Fe/H] at that redshift, according 
to the derived relations. In addition, in each galaxy the GC 
metallicity is assumed to have a gaussian distribution around 
the mean value, with a ~ 0.3 dex. 

4. RESULTS 

Fig. <J2J shows the globular cluster metallicity distribution, 
in a galaxy of mass Mo = 10 11 M© at redshift z = 0, with pa- 
rameters Mi /M = 0.3 and M SF /M = (i.e. 70% of the final 
stellar mass come from accreted satellites, and there is no ad- 
ditional star formation in the main progenitor after the GC 
formation), averaged over N = 1000 Monte Carlo realisations 
of the galaxy formation history. The thick black line shows 
the total distribution, while the red line shows the distribution 
for the clusters that were formed locally in the main progeni- 
tor at z = Zj n ~ 2, and the blue line shows the distrubution for 
the clusters formed in satellites at an epoch z ~ 3 -4, and that 
merged with the main progenitor. The dashed lines + shaded 
areas show the values of the metal-rich and metal-poor GC 
metallicity [Fe/H] and their la uncertainties for a galaxy of 
mass Mo from the Peng et al. (2006) relations. 

The model galaxy shows a sharp bimodality in the globu- 
lar cluster metallicity distribution. The metal-rich peak of 
the metallicity distibution is entirely dominated by local GCs, 
formed in the main progenitor at z ~ 2. The metal-poor peak 
is entirely dominated by satellite GCs, accreted via the hierar- 



chical assembly. The positions of both peaks are consistent 
with the observed galaxy stellar mass - GC metallicity re- 
lations by Peng et al. (2006) for metal-rich and metal-poor 
globular clusters. 

The number of GCs in various realisations of this galaxy 
scatters around the value Tn ~ 6 interpolated from Peng et 
al. (2008) for a galaxy of mass M star ~ 10 n M Q , staying in 
the (rather large, ±5) observed scatter limits. For a given set 
of history parameters Mi /Mo and Msf/Mo, the final value 
of Tn depends mainly on the assumption about the minimum 
stellar mass of a galaxy that can form globular clusters (here 
Miimit = 1O 9 M ). Note that a variation of a factor 10 in the 
mass limit, such that M]; m i t = 10 8 M Q , yields a value Tn ~ 50, 
one order of magnitude off the Peng et al. (2008) relation. 

The positions of the peaks in the GC metallicity distribution 
are determined by the galaxy [Fe /H] - M star relation, given 
that we can constrain the ages of the GCs from observations, 
and under the assumption that the GC metallicity is connected 
to the instantaneous metallicity of the galaxy where they were 
formed. But is the bimodality just a result of these choices, 
or is it an intrinsic feature of our mass assembly scenario? 
Fig. © shows the GC metallicity distribution of a galaxy 
characterised by (M = 10 11 M Q , Mi/M = 0.5, M SF = 0) , this 
time under a very conservative assumption: the GC metallic- 
ity in all the objects in the merger tree (main progenitor and 
satellites) is just taken from the average [Fe/H]Gc~M star re- 
lation of Peng et al. (2006) (thick black line in FigdJ. 

Although both the local and the accreted GCs obey the 
same average relation, they are still separated in metallicity, 
the distribution of which shows two distinct peaks, albeit at 
the wrong values. The reason why the metallicity bimodal- 
ity is still present is that the hierarchical mass assembly is 
governed by the halo merger rate, which greatly suppresses 
merger events of high mass ratios (larger than 1:10) (as ev- 
ident in Fakhouri et al. 2010), so that it is highly improb- 
able that a galaxy merges with objects of similar mass, and 
therefore similar metallicity. This feature alone is what drives 
the bimodality in the GC metallicity distribution. Therefore, 
a metallicity bimodality in the GC population is a direct 
prediction of the hierarchical clustering scenario. 

Notice also that the metal-poor peak in Fig. (© is almost at 
the right value of [Fe/H], while the metal-rich peak is off by 
~ 0.5 dex towards the metal-poor side. The slope of the aver- 
age [Fe/H]cc - M stal - relation in Fig Q] suggests that the num- 
ber of metal-poor GCs is highly dominant in low-mass galax- 
ies. This happens because these are intrinsically metal-poor 
galaxies; in addition note that, as their stellar mass is small, 
in their assembly history they are only able to accrete smaller 
objects that are devoid of globular clusters (given the exis- 
tence of Miimit), therefore their GC population is not bimodal, 
and their average metallicity peaks exactly where the metal- 
poor peak is located. On the other hand, the more massive 
a galaxy is, the richest its assembly history is, with a merger 
tree with enough mass range to sustain a varied secondary GC 
population, so its GC population is more likely to be bimodal. 
Therefore, a massive galaxy always has a secondary, metal- 
poor GC population, and the average GC metallicity deviates 
from both peaks. This point is addressed in the next Figure. 

Fig. © illustrates the difference in the GC total metal- 
licity distribution of galaxies with final stellar masses Mo = 
10 n ,5 x 10 1() , 10 1() ,5 x 10 9 M Q , all characterised by the his- 
tory parameters Mi /Mo = 0.3, M S f/M = 0. The bimodality 
in the GC metallicity distribution is evident at all masses in 
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FIG. 2. — The GC metallicity distribution of a galaxy of stellar mass Mo = 10" Mq. The history of this galaxy is characterised by the parameters M| /Mo = 0.3, 
Msf/Mo = 0. The local GC metallicity is sampled from a gaussian distribution centered around the galaxy [Fe/H] — M s tar relation at redshift z ~ 2, while the 
metallicity of GCs accreted from satellites is centered around the galaxy [Fe/H] — Mstar relation at redshift z ~ 3.5. Black line: total GC metallicity distribution; 
red line: metallicity of local GCs; blue line: metallicity of GCs accreted from satellites. Dotted lines and shaded areas: values of [Fe/H] for the metal-rich and 
metal-poor GC populations of a galaxy of M stiU = Mo from the relations of Peng et al. (2006) and corresponding scatter. 




FIG. 3. — The GC metallicity distribution of a galaxy characterised by 
(M = lei 1 M©, M! /Mo = 0.5, M S f = 0), where the GC metallicity of both 
the local and the accreted GC populations are sampled from a gaussian distri- 
bution centered around the average globular cluster [Fe/HJoc ~ Mstar relation 
from Peng et al. (2006: black line of Fig. 0. The dotted lines and shaded 
areas represent the values of [Fe/H] for the metal-rich and metal-poor GC 
populations of a galaxy of M stal = Mo from the relations of Peng et al. (2006) 
and corresponding scatter. 
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FIG. 4. — The total GC metallicity distribution of galaxies of varying final 
stellar mass, all characterised by the parameters Mi /Mq = 0.3, Msf /Mq = 0. 



this mass range. As expected, the more massive a galaxy is, 
the richest is its GC population, in both the metal-rich and the 
metal-poor component. However, notice that, although the 
history parameters are the same in all cases, the relative con- 
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tribution of the two peaks varies, with the metal-poor peak be- 
coming less and less significant relative to the metal-rich peak 
for lower-mass galaxies, in accord with Peng et al. (2008) 
and Shapiro et al. (2010). In the lowest mass bin, the rela- 
tive height of the metal-rich and metal-poor peaks is reversed; 
given that the mass limit for GC formation is 10 9 M Q , this 
galaxy is for the most part accreting satellites that don't con- 
tribute to the GC population, with the rare exception of major 
mergers (in this case, M sat > 10 9 M Q ). If we assume that the 
globular cluster formation is hampered in low-mass galaxies, 
i.e. that galaxies below a mass threshold cannot produce glob- 
ular clusters, then the model predicts that the GC metallicity 
bimodality ceases to exist slightly above that mass threshold. 
In such galaxies, the GC population is unimodal and entirely 
composed of locally-formed GCs. On the other hand, the 
overall GC metallicity decreases following the galaxy mass, 
and as a result, the GC population in low mass galaxies is 
metal-poor, again in accord with Peng et al. (2008). 

The galaxy assembly history determines the fraction of the 
final mass that is accreted from the merger tree, and there- 
fore the fraction of globular clusters that are formed outside 
the main progenitor and which we have shown to compose 
the metal-poor peak. Fig. (|5]l shows the relative height of the 
metal-rich and metal -poor peaks generated in different assem- 
bly histories, parameterised as Mi /Mo. For a galaxy of final 
stellar mass Mq = 10 11 M , the panels from left to right show 
the GC metallicity distribution for Mi /M = (0.2,0.3,0.5,0.6) 
respectively. In all cases, Msf = 0. As expected, a galaxy with 
a poor merger history (such as the case Mi /Mo = 0.6 for in- 
stance) shows a GC metallicity distribution dominated by the 
local metal-rich population. The model therefore predicts that 
the presence of a very strong metal-rich GC component is a 
sign of a sparse merger history. For a given galaxy mass, the 
richness of the merger tree depends on environment; hence 
the model predicts that galaxies in low-density environments 
have, for a given mass, a GC population that is more metal- 
rich dominated than galaxies of the same mass living in the 
centre of clusters. 

So far we have analysed the simplified case of galaxies with 
Msf = 0. However, for the majority of galaxies the star forma- 
tion does not stop at z ~ 2, and a significant part of the final 
stellar mass is formed at later times. In this case, a signifi- 
cant fraction of the galaxy stellar mass is not associated with 
formation or accretion of globular clusters. To account for 
this stellar component, we vary the value of the star formation 
history parameter Mg F /Mo. Fig. (|6]l shows the GC metallicity 
distribution for a galaxy of final mass Mo = 10 11 M , where 
the stellar mass is contributed by 1) the main progenitor at 
the epoch of GC formation in proportion of Mi /Mo = 0.3 (the 
local GC population is associated with this component), 2) 
stars formed locally in the galaxy after the epoch of GC for- 
mation, in quantity Msf/Mo = (0.2,0.3,0.5,0.6) {panels from 
left to right), and 3) stars accreted from satellites, in quantity 
M -Mi -M S f = M sat . Note that a higher value of M S f/M 
implies a smaller value of M sat , i.e. a poorer merger history. 

Fig. © shows that, as the galaxy growth becomes more 
dominated by local star formation and the contribution of the 
mass accreted by satellites is smaller, the GC population be- 
comes more and more dominated by local globular clusters, 
even if most of the stellar component is not directly associ- 
ated with the globular clusters themselves. The model pre- 
dicts that, in galaxies with an active star formation history 
after the GC formation (i.e. at z < 2), the relative strenght 



of the metal-rich and metal-poor peaks of the GC metallicity 
distribution is biased towards the metal-rich GCs, for a given 
galaxy mass. If we consider the star formation history as as- 
sociated with morphology, then the model predicts that, for 
a given galaxy mass M and total number of GCs, late-type 
galaxies have a GC metallicity distribution with a stronger 
metal-rich peak than early-type galaxies. 

Note that the results in Figs. ([5] |6]l show that, given the 
final mass of the galaxy Mo, the final number of globular 
clusters in the galaxy depends on the value of the assembly 
history parameter Mi /Mo and the star formation history pa- 
rameter Msf/Mo. The large scatter in the value of the GC 
frequency per unit mass Tn for a given galaxy mass seen in 
Peng et al. (2008) is likely to be due to the variety of his- 
tories for galaxies in each mass bin. The final value of Tn 
decreases for a decreasing value of M sat . Note that the as- 
sembly and star formation history parameters have instead no 
effect on the position of the peaks, which are entirely deter- 
mined by the evolution of the galaxy [Fe/H]-M star relation. 
The number of GCs, together with the relative abundance of 
the metal-rich and metal-poor components, can therefore be 
used to constrain the assembly and star formation history of 
the galaxy. 

As discussed in the Introduction, the scenario in which 
globular cluster are generally formed in gas-rich mergers can- 
not reproduce the properties and the scaling relations of the 
GC population. However, if globular clusters indeed form 
during violent bursts of star formation, it is physically possi- 
ble that some of them indeed are formed in mergers at all red- 
shifts, a fact that would explain the presence of intermediate- 
age or young GCs in some galaxies (Kissler-Patig et al. 1998, 
Puzia et al. 2005, Strader et al. 2003, 2004b, Woodley et 
al. 2010, Brodie & Strader 2006 and references therein). 
Fig. (Q shows the effect of a gas-rich merger event, where 
new GCs are formed, on the GC metallicity distribution. In 
this example, in the same galaxy portrayed in Fig. (O, we 
introduce a gas-rich merger event that triggers the formation 
new GCs, in quantity ~ 30% of the local GC population of 
the main progenitor, with intermediate metallicities peaked 
around [Fe/H] ~ -0.8 with a gaussian of width a = 0.2 dex 
{green line). This plot shows that the creation of new GCs 
in merger events introduces a stochastic variation of the GC 
metallicity distribution, that leads to the formation of ter- 
tiary peaks, in positions determined by the metallicity of the 
gas perturbed/carried by the merger. The number of newly- 
formed GCs depends on the available gas mass and the star 
formation rate in the merger-triggered bursts, as well as the 
efficiency of GC formation vs star formation. 

It is clear from this plot that, if we consider the formation of 
GC in gas-rich merger events, the GC metallicity distribution 
becomes more complex. The stochasticity of such events al- 
lows for any shape of the final metallicity distribution: a pro- 
longued history of gas-rich mergers contributes to the dilution 
of the bimodality. Such mechanism can explain the number of 
'exotic' GC metallicities distributions found by a number of 
authors, including Blom et al. (2012) and Peng et al. (2006), 
with a number of galaxies that either show one or more ter- 
tiary peaks, or a non-bimodal GC metallicity distribution (it 
should be noted that such a scenario needs to be confirmed 
with dynamical analysis; for instance, Blom et al. 2012 show 
data of a galaxy with an intermediate-metallicity GC subpop- 
ulation that rotates with the main body of the galaxy). Note 
that major gas-rich mergers are good candidates to provide 
very intense bursts of star formation, during which new globu- 
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FIG. 5. — The variation of the GC metallicity distribution for different assembly histories, parameterised by Mi/Mo. From left to right: Mi /Mo = 
(0.2,0.3,0.5,0.6). In all cases, the final galaxy stellar mass is Mo = 10" Mq, and Msp = 0. The dotted lines and shaded areas represent the values of 
[Fe/H] for the metal-rich and metal-poor GC populations of a galaxy of M sta i = Mq from the relations of Peng et al. (2006) and corresponding scatter. 
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FIG. 6. — The variation of the GC metallicity distribution for different star formation histories, parameterised by Msf/Mo. From left to right: Msf/Mo = 
(0.2,0.3,0.5,0.6). In all cases, the final galaxy stellar mass is Mo = 10" Mq, and Mi /Mo = 0.3. The dotted lines and shaded areas represent the values of 
[Fe /H] for the metal-rich and metal-poor GC populations of a galaxy of M slal = Mo from the relations of Peng et al. (2006) and corresponding scatter. 
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FIG. 7. — The GC metallicity distribution of a galaxy of mass Mo = 
10 11 M with Mi /M = 0.3 ( see Fig. (3))> with the inclusion of a major gas- 
rich merger event in the assembly history. The merger triggers the formation 
of a tertiary population of GCs, in number N = 0.3 X N(Mi) (i.e. 30% of the 
main progenitor's local GCs) and intermediate metallicities peaked around 
[Fe/H] ~ -0.8 (green line). 



lar clusters can be formed. If a galaxy undergoes an assembly 
history devoid of any gas-rich mergers, it is hard to envisage 
another mechanism that is able to provide a strong enough 



perturbation of the gas in the galaxy to trigger very intense 
bursts of star formation (as per Shapiro et al. 2010), which 
can create a tertiary GC population. Therefore, we can con- 
sider tertiary peaks in the GC metallicity distribution as clear 
signatures of major gas-rich merger events in the past history 
of the galaxy. 

The total number of GCs and the relative height of the 
metal-rich and metal-poor peaks depend on the galaxy mass 
and the galaxy assembly and star formation history, while 
the positions of the peaks depend on the determination of the 
galaxy mass-metallicity relation as a function of redshift. To 
explore this point further, we study the GC metallicity distri- 
bution resulting from the [Fe /H] - M stal - relation obtained in 
Section 3, under different prescriptions. In particular, we de- 
scribe two examples of variations of the recipes described in 
Section 3, that affect 1) the normalisation and 2) the evolution 
of the galaxy [Fe/H] - M star relation: 1) is a variation of the 
value for the solar 12 + log(0/H), and 2) is a variation of the 
prescription for the evolution of the galaxy [a/Fe] value. 

Fig. ^ shows the GC metallicity distribution of a galaxy 
characterised by (Mo = 10 n M Q , Mi/M = 0.3, M SF = 0) 
(right panels) resulting from two different sets of galaxy 
[Fe/H] - M star relations as a function of redshift (left panels), 
which are derived following the model of Section 3. In the 
upper left panel, the galaxy [Fe/H] — M star relation is obtained 
from the evolution of 12 + log(0/H), but assuming that the 
galaxy [O/Fe] does not evolve with redshift, but has a con- 
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FIG. 8. — The dependence of the GC metallicity distribution on the evolution of the galaxy mass-metallicity relation. Upper left panel: the galaxy [Fe/H]-M slar 
relation as a function of redshift (as per Section 3, compare with Fig. (T}), in the case of no [O/Fe] evolution. Upper right panel: the corresponding GC 
metallicity distribution. Lower left panel: the galaxy [Fe/H]— M sta r relation with the [O/Fe] evolution described in Section 3, but with the oxygen solar value set 
as 12 + log(0/H) = 8.9. Lower left panel: the corresponding GC metallicity distribution. In the right-hand side panels, the portrayed galaxy is characterised by 
(Mo = 10' 1 Mq , Mi /Mo = 0.3, Msf = 0); the dotted lines and shaded areas represent the values of [Fe/H] for the metal-rich and metal-poor GC populations of 
a galaxy of M slar = Mq from the relations of Peng et al. (2006) and corresponding scatter. 



stant value of ^0.1 like in the local universe (with a solar 
oxygen value 12+log(0/H) = 8.66 as in our fiducial relation). 
In the upper right panel, the resulting galaxy GC metallicity 
distribution shows a somewhat diminished bimodality, and the 
positions of the peaks are definitely off the observed values 
obtained by Peng et al. (2006; dotted lines and shaded areas). 
Both peaks are centered around too high metallicities, and the 
problem is worse for the metal-poor peak; if we assume that 
the galaxy [O/Fe] does not evolve with redshift but mantains 
the local value, then we are overestimating the galaxy metal 
content by a factor that is proportional with redshift, and the 
globular clusters formed in small galaxies at high redshifts are 
the ones that are affected the most. 

The lower left panel of 'Fig. dS) portrays the galaxy [Fe/H]- 
M star relations as obtained in Section 3, with the same evolu- 
tion of [O/Fe] used so far in this work, but with a different 
solar oxygen value, 12 + log(0/H) = 8.9 (Liu et al. 2008). In 
this case, the evolution of the [Fe /H] — M sta , relation is not af- 
fected, but its normalisation has changed. The resulting GC 



metallicity distribution does not change in shape (the relative 
height and position of the peaks is the same) but both peaks 
are shifted towards lower [Fe /H] values, because the oxygen 
content [O/H] yields a lower total [Fe/H] content. 

The model predicts that the positions of the metal-rich and 
metal-poor peaks of the GC metallicity distribution are exclu- 
sively dependent on the galaxy mass-metallicity relation as a 
function of redshift. In the case where globular clusters are in- 
deed fossil records of the metallicity of their parent galaxy at 
the time when they formed, then through this model the GC 
metallicity distribution, and in particular the position of the 
metal-rich and metal -poor peaks, can be used to constraint the 
evolution of the galaxy [Fe /H] - M star relation and the evolu- 
tion of the galaxy [O/Fe]. 

5. DISCUSSION 

The results described in the previous Section show that the 
hierarchical galaxy assembly directly predicts a metallicity bi- 
modality of the globular cluster populations in galaxies. Pre- 
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vious studies like Shapiro et al. (2010) on the contrary, argued 
that the hierarchical mass assembly would blur the correla- 
tions between galaxies and their GCs. The model presented 
here predicts globular clusters to have a bimodal metallicity 
distribution, with the metal-rich and the metal-poor peaks fol- 
lowing the [Fe/H]Gc - M star relations described in Peng et al. 
(2006), and determined by the evolution of the galaxy mass- 
metallicity relation. The model also predicts that the relative 
strenght of the metal-rich and metal-poor peaks depends on 
the assembly and star formation history of the galaxy. More- 
over, the model predicts that the metallicity bimodality dis- 
appears for masses below M star ~ 1O 9 M and for redshifts 
z > 2. 

The main mechanism at work in producing the GC bi- 
modality is the existence of the galaxy mass-metallicity rela- 
tion, coupled with the natural behaviour of hierarchical mass 
assembly to strongly favour minor mergers over major merg- 
ers. In this scenario, all globular clusters form with the same 
mechanism and share the metallicity of their parent galaxy; 
the hierarchical build-up of galaxies then assembles the GC 
populations so that the GCs of a satellite become part of 
the accreted, metal-poor GC component of a bigger galaxy 
of higher metallicity (in a generalisation of the scenario of 
Shapiro et al. 2010). Note that this mechanism still works 
in the case of satellites with bimodal GC distributions: both 
the metal-rich and the metal-poor subpopulations of a satel- 
lite are metal-poor compared to the central galaxy, and will 
contribute to its final metal-poor GC population. 

Under the assumption that GC formation is a rapid process, 
closely linked to the properties of the galaxy at the time of the 
event, the fact that GCs of different metallicities form in dif- 
ferent galaxies, at slightly different times, is enough for the hi- 
erarchical galaxy formation scenario to naturally produce the 
correct GC metallicity distribution. The observational facts 
imply this very clearly. For instance, a comparison of the 
metal-poor globular cluster [Fe /H]gc _ M star relation with the 
galaxy mass-metallicity relation (blue straight line in Fig. (T|)) 
shows that, for any given galaxy mass above M star > 10 9 M Q , 
the galaxy is always more metal-rich that its own metal-poor 
globular clusters (unless they were formed at z > 4). This im- 
plies that these globular clusters must have formed in smaller 
galaxies and they were then accreted. Note that the alternative 
scenario of a multi-phase GC formation in each single galaxy 
(like for instance Beasley et al. 2002), where metal-poor GCs 
are formed first and metal-rich GCs are formed later in merg- 
ers, necessarily implies some form of systematic segregations 
of metals and an ad hoc shut-down of the formation of the 
metal-poor component, in order to produce both the correct 
metallicity and the bimodality. In fact, if a galaxy experiences 
a prolongued phase of local GC formation, the bimodality is 
destroyed by the galaxy chemical self-enrichment. 

The results presented in this work are based on theoreti- 
cal merger trees extracted from the Millennium simulation. 
Among the cosmological parameters, the value of a% can af- 
fect these results, in the sense that a lower clustering power 
would produce a delayed mass accretion and sparser merger 
histories. While the bimodality would remain unchanged, the 
height of the metal-poor peak would be affected. However, 
the magnitude of the effect would be smaller than the scatter 
in the observed values of TJv and the scatter between differ- 
ent Monte Carlo runs, and the results presented in this work 
would remain unaffected. In any case, this would be an inter- 
esting avenue of investigation. 



The bimodality in the GC metallicity distribution is solid 
against variations of our initial conditions, such as any as- 
sumption about the metallicity we assume for the galaxy or 
the GC themselves (Figs. 12 [3] EJ. Thorugh hierarchical 
galaxy assembly, it is rather hard to produce a GC metallicity 
distribution that is not bimodal. Nonetheless, the position of 
the peaks is not recovered correctly if we assume that the local 
and accreted GC populations form at the same time, given the 
galaxy mass-metallicity relations described in this work. This 
is in accord with observations that determine the metal-poor 
GCs to be 1-2 Gyr older than the metal-rich GCs. A second- 
order sophistication of this model would be to introduce an 
analytic relation between the galaxy mass and the epoch GC 
formation, but in reality this epoch is likely to vary with envi- 
ronment and the fluctuations of the local star formation rate. 
A scatter in the epoch of formation would mirror a scatter in 
the galaxy metallicity, via the evolution of the galaxy mass- 
metallicity relation; this effect has been mimicked in this work 
by the introduction of a scatter in the GC metallicity of the 
satellites. 

In this model, old globular clusters are considered fossil 
records of the galaxy where they were formed, and their fi- 
nally metallicity distribution is exclusively a result of the hi- 
erarchical galaxy assembly. The factors that affect the final 
metallicity distributions are 1) the evolution of the galaxy 
mass-metallicity relation, which completely determines the 
position of the metal-rich and metal-poor peaks, and 2) the 
merger history (therefore the environmental density) and the 
star formation history, which completely determine the rela- 
tive strenghts of the metal-poor and metal-rich peaks. These 
constraints are independent and orthogonal, with no degener- 
acy between them. 

Provided we know the GC ages, we can use this model to 
test and constrain the evolution of the galaxy [Fe /H] - M star 
relation and the evolution of [a/Fe], through the positions 
of the metal-rich and metal-poor peaks. In order to repro- 
duce the Peng et al. (2006) [Fe/H] G c -M stal - relations, the 
model favours a value of [a/Fe] = 0.5 at redshift z = 4 and 
a linear evolution down to [a/Fe] =0.1 at z = 0, a solar 
value 12 + log(0/H) ~ 8.66, and the evolution of the galaxy 
[Fe/H] -M stal relation plotted in Fig. (HJ. Note however that, 
even if we can rely on spectroscopy for the determination of 
the globular cluster metallicity with good precision, there are 
still significant uncertainties on the GC age determination. In 
this work we have used the observed estimates for the aver- 
age GC ages, but if we consider their uncertainties, combined 
with the current size of the uncertainly on the galaxy mass- 
metallicity relation, then there is a substantial degeneracy be- 
tween redshift and galaxy metallicity in the determination of 
the positions of the metal-rich and metal-poor peaks, an issue 
that will be solved with higher precision observations from 
the next-generation instruments and surveys. 

Unrelated to the particular values of the metal-rich and 
metal-poor peaks, the model predicts the relative height of the 
peaks to give an insight into the assembly and star formation 
history of the galaxy: for a given galaxy mass, a dominant 
metal-rich peak indicates a quiet merger history and/or a sig- 
nificant growth of the galaxy through local star formation, as 
opposed to a mass growth driven by accretion of satellites. 
Therefore, the relative abundance of metal-poor and metal- 
rich GCs is correlated with environmental density and mor- 
phology, with isolated late-type galaxies of given mass Mo 
being more metal-rich GC dominated than early types of the 
same mass living in dense environments. The current state of 
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observations of GC populations in late-type galaxies does not 
yet provide definitive constraints, but this is certainly an area 
where more numerous and improved observations are called 
for. This is of particular interest when attempting to con- 
strain the galaxy star formation history, in that this method 
is complementary to the SED-fitting technique, which is af- 
fected by a number of systematic uncertainties (Tonini et al. 
2009, 2010, 2011, 2012, Henriques et al. 2011, Maraston et 
al. 2010, Pforretal. 2012). 

The hierarchical galaxy mass assembly naturally leads to 
bimodality in the globular cluster metallicity distribution. A 
bimodal distribution is sign of a two-phase galaxy formation, 
with an intense dissipative phase that leads to the formation of 
the core of the galaxy and the local GCs, and a second phase 
of accretion of satellites and the acquisition of a secondary GC 
population. This regularity breaks down at very small masses, 
when galaxies are not massive enough to accrete GCs from 
their merger tree: their satellites are in fact so small, that they 
cannot form their own GCs. In this case the GC distribution 
is unimodal. On the other hand, the smaller a galaxy is, the 
least probable it is that it actually has an extended merger tree, 
so even when the few satellites contribute with globular clus- 
ters, the metal-poor peak is subdominant or negligible. Note 
that this behaviour is not in disagreement with the Tn - M star 
relation presented in Peng et al. (2008). For masses below 
10 9 M©, even 1 globular cluster will yield a value Tn ~ 10. 
This shows that our limit mass for the production of globular 
clusters, Mii m ; t ~ 10 9 , is a realistic prediction. 

In relation to the number of globular clusters per galaxy 
mass Tn predicted by the model, the values we obtain are in 
the range of the scatter observed by Peng et al. (2008). We 
argue that in observations such scatter arises from the variety 
of assembly and star formation histories that generated the 
galaxies in the sample. In the model, the scatter between dif- 
ferent Monte Carlo runs is of the same order of magnitude, 
and we argue that such a scatter generates from the scatter 
in the mass function of the merger tree, incorporating all the 
possible assembly paths to build up the final galaxy mass. 
The final value of T N depends on the assembly and history 
parameters Mi /Mo and M SF /M , where richer merger histo- 
ries strenghten the metal-poor peak and raise the total number 
of globular clusters. 

The model also predicts that the GC bimodality is a function 
of redshift. The more time a galaxy has to accrete satellites 
after the GC formation, the more rich its secondary population 
will be. Therefore, we expect the GC metallicity bimodality 
to disappear by redshift z ~ 2 and above. 

6. SUMMARY AND CONCLUSIONS 

We have performed a series of Monte-Carlo simulations of 
the assembly history of galaxies and studied the metallicity 
distribution of their globular cluster systems. To do so, we 
built a model for the assembly of the globular cluster popula- 
tion following the hierarchical galaxy assembly, based on em- 
pirical scaling relations such as the galaxy mass-metallicity 
relation [0/H]-M star , and on the observed galaxy stellar mass 
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function up to redshift z ~ 4. We also made use of the theoret- 
ical merger rates as a function of mass and redshift from the 
Millennium simulation, to build merger trees for a set of final 
galaxies. 

By determining the galaxy [Fe/H]-M star relation for all 
galaxies in each merger tree, and by assuming that globular 
clusters share the metallicity of their original parent galaxy 
at the time of their formation, we populated the merger tree 
with globular clusters. The hierarchical assembly of the final 
galaxy creates a globular cluster population composed by the 
local GCs formed in the main progenitor and those accreted 
from the merger tree. We conclude that: 

• the final globular cluster metallicity distribution is in gen- 
eral bimodal; the GC metallicity bimodality is a direct predic- 
tion of the hierarchical clustering scenario; 

• the metal-rich peak of the GC metallicity distribution is 
composed of globular clusters locally formed in the main pro- 
genitor, while the metal-poor peak is composed of globular 
clusters accreted from the satellites that compose the merger 
tree. At all times GCs in satellites are more metal-poor than 
GCs formed in the main progenitor due to the existence of the 
galaxy mass-metallicity relation; both the metal-rich and the 
metal-poor subpopulations of a satellite will contribute to the 
metal-poor GC population of the main galaxy; 

• the positions of the metal-rich and metal-poor peak 
depend exclusively on the evolution of the galaxy mass- 
metallicity relation [Fe/H] -M star ; we are able to constrain 
such evolution and predict that the galaxy [O/Fe] evolves lin- 
early with redshift from a value of ~ 0.5 at redshift z ~ 4 to a 
value of ~ 0. 1 at z = 0; 

• the relative strenght of the metal -rich and metal -poor peak 
depends on the assembly and star formation history of the 
galaxy. The model predicts that, for a given galaxy mass, 
galaxies with a poor merger history, such as galaxies form- 
ing in low density environments, and/or galaxies with a pro- 
longued star formation history (after the epoch of GC forma- 
tion) that contributes most of the galaxy mass, such as late- 
type galaxies, will have a globular cluster population domi- 
nated by the metal-rich component. On the other hand, galax- 
ies of the same mass but with an intense merger history, such 
as early-type galaxies and/or galaxies living in dense environ- 
ments, will have a globular cluster population with a larger 
metal-poor component; 

• the model predicts that the globular cluster metallicity bi- 
modality disappears at galaxy masses around 10 9 M©; more- 
over, the model predicts that the bimodality is progressively 
less pronounced at higher redshift, and disappears around red- 
shift z - 2. 
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